A Non-Disordered Glassy Model with a Tunable Interaction Range 
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We introduce a non-disordered lattice spin model, based on the principle of minimizing spin-spin correlations 
up to a (tunable) distance R. The model can be defined in any spatial dimension D, but already for D = 1 and 
small values of R (e.g. R = 5) the model shows the properties of a glassy system: deep and well separated 
energy minima, very slow relaxation dynamics, aging and non-trivial fluctuation-dissipation ratio. 

PACS numbers: 75.50.Lk,05.50+q,75.40Mg 



In a supercooled liquid the viscosity increases abruptly by 
several order of magnitude in a narrow temperature range, and 
eventually undergoes a dynamical arrest, that can be observed 
on any accessible time scale: this is the essence of a fasci- 
nating physical phenomenon, the so called "dynamical glass 
transition". A vast scientific literature has been dedicated to 
its study: see for interesting reviews on the subject. A 
theoretical understanding of this effect must be based on a re- 
liable modeling of the underlying material: ideally one would 
like to have at hand very simple models that reproduce the 
main features of glass-former liquids. From the point of view 
of a numerical approach, off-lattice simulations are extremely 
costly in term of computational resources: because of that lat- 
tice models where each node of the interacting network con- 
tains some degrees of freedom (a binary spin variable in the 
simplest case) play an important role. 

An analytic approach needs some reasonable approxima- 
tion: mode coupling theory |ll|,|2|], for example, helps to shed 
some light on the problem. At the mean-field level, and more 
precisely on a fully connected lattice, the solution of the disor- 
dered p-spin model 03[] is now very well understood [4], and it 
turns out to be equivalent to the mode-coupling theory (based 
on systems where the Hamiltonian does not include quenched 
disorder). The main prediction of these mean-field theories is 
that below the dynamical glass transition T g , which is higher 
than the thermodynamical critical point Tk, the relaxation dy- 
namics is not able to bring the system to equilibrium in any 
sub-exponential time (in the system size). Consequently the 
system relaxes to the so-called threshold energy, which is 
higher than the equilibrium energy, and the off-equilibrium 
dynamics shows aging on any measurable time scale lllql . 
During the aging dynamics the fluctuation-dissipation ratio is 
different from the one expected in equilibrium, and its behav- 
ior presents peculiar and distinctive features. 

The mean field scenario is well characterized and well un- 
derstood, but it is still unclear how to adapt it to real systems. 
In finite dimensional systems the lifetime of metastable states 
is limited: eventually, during the aging dynamics, a bubble 
of the equilibrium state will nucleate and will grow up to the 
system size. Nonetheless the time for nucleating and grow- 
ing the equilibrium phase may be extremely large, especially 
close to a critical point. Moreover the Hamiltonian of a real 



glass-former does not contain any quenched disorder: this is 
a crucial difference from the p-spin model. The frustration, 
which is the main ingredient for the slow relaxation, is self- 
induced by the relaxation dynamics: a realistic model of a 
glassy system should contain no quenched disorder. 

The "kinematic models" [5] are glassy models with no 
quenched disorder, where the evolution is governed by a spe- 
cific dynamical rule, but it does not correspond to a relaxation 
on a well-defined energy landscape: they are interesting but 
they cannot undergo a true thermodynamical phase transition, 
and we will not consider them in the following. 

Not so many "Hamiltonian glassy models" without quen- 
ched disorder are available. Shore, Holzer and Sethna ||6(] in- 
troduced and analyzed a model with a competition in the inter- 
action, and a tiling 2D model. Newman and Moore |7|] have 
discussed a 2D model with a triangle 3-spin interaction, that 
can be solved exactly: it does not show a sharp glass transition 
but undergoes a severe slowing down. Biroli and Mezard [8] 
defined a very simple lattice glassy model, where the legal 
positions of the particles are restricted by hard "density con- 
straints": the model is versatile, since it does not depend on 
the detailed feature of the underlying lattice, but its energy 
landscape is somehow drastic, in the sense that a configura- 
tion is either legal (E = 0) or forbidden (E — oo). Cavagna, 
Giardina and Grigera [9] have discussed a 2D model based 
again on competing interactions (respectively with 4 and 5 
spins). It is clear that enlarging this collection would be ap- 
propriate: some of these interesting models are indeed strictly 
two dimensional or depend on the detailed lattice structure. In 
some other cases one can observe a very slow domain growth, 
but once the time is appropriately rescaled the growth process 
does not differ qualitatively from the dynamics of the pure 
Ising model. 

In this note we define and analyze a new non-disordered 
glassy model, that can be defined on any lattice structure, in 
any dimension, and is based on a simple physical principle: 
the minimization of correlations. Our model is also, following 
the route of its mean field predecessor, a good candidate for 
providing coding for an effective and secure communication. 

We are inspired from the Bernasconi mean field model Jiol 
1 111 , where one is interested in finding the assignment to N 
Ising spins defined on a linear lattice that minimizes the sum 
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= j^j EdLi 1 (T,f=o d 1 a 3 a j+d) of the squared spin- 
spin correlations. We have used here open boundary condi- 
tions, but the Bernasconi model, like our model, is also inter- 
esting when defined with periodic boundary conditions. The 
Bernasconi model has a very rough energy landscape |l2 , 131. 
with deep minima separated by extensive energy barriers, 
making the search for global minima (low auto-correlation bi- 
nary sequences) a very difficult task. The theoretical analy- 
sis of this model predicts a thermodynamical phase transition 
with one step of replica symmetry breaking (1RSB), in the 
same universality class of the p-spin model, preceded by a 
dynamical glass transition. Extensive numerical simulations 
have shown that the energy relaxation stops before reaching 
the ground state (GS) energy, and the aging regime persist for 
extremely long times. For these reasons the Bernasconi model 
is the perfect prototype of a glassy mean field model. 

We define our model by adopting the principle of minimiz- 
ing spin-spin correlation functions, but using an interaction 
that is local in space: in D dimensions the Hamiltonian reads 



max d in IZ(x) j couples in IZ(x) 
d— 1 \ at distance d 



(l) 



where A C IP is a finite volume of cardinality TV, TLix) 
is a hypercube of size RP (or its intersection with A if 
open boundary conditions are used) centered around site x, 
and Nn.b. is a normalization constant that guarantees good 
R . —y oo and N — > oo limits. The sum over d is for dis- 
tances going from one up to the maximum distance contained 
in lZ(x). One can also define the model by only considering 
correlations on the D axis: from the point of view of compu- 
tational cost this is a far better choice. 

With eq. (QJ we are aiming at minimizing correlations in 
blocks of linear size R. Since correlations at short distances 
are typically the strongest, the overall effect is to have low- 
energy configurations showing very weak correlations on all 
length scales. This seems to us a very solid first principle to 
build a glassy model: glass-formers show no long range order 
in the two-point correlation functions and we are somehow 
enforcing this condition in the Hamiltonian. The tunable in- 
teraction range R is a novel and very useful feature of our 
model. The R — > oo limit gives back the Bernasconi model 
in D = 1, and for D > 1 provides new and potentially inter- 
esting mean-field models. As for the Bernasconi model here 
we can get a good basis towards effective coding: the intro- 
duction of new, free parameters (D and R) that are unknown 
to the observer could be of further help. 

Let us look in better detail to the D = 1 version of the 
model with open boundary conditions. The Hamiltonian reads 
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Y,C{d,i,R)\ (2) 



where C(d,i,R) = £ 



i+R-1- 



' aj<jj+d, and R is the tunable 



interaction range. We will show, for example, that the R = 5 



model presents to a very large extent the glassy phenomenol- 
ogy we have described above (as opposed, for example, to the 
R = 6 model). This makes clear that already the D — 1 case, 
where a thermodynamical transition is forbidden, is very in- 
teresting for the physics of structural glasses. 

Finding and analyzing low energy configurations is our first 
task. For R = N (the Bernasconi model) finding a GS is 
difficult in practice and requires a time growing exponentially 
with N. For finite R the model can be solved in a time of order 
0(N2 2 ( R ~ 3S> ) by transfer matrix methods, using the variables 

T s ; = (Tj_l(Tj + l. 

We have computed GS and first excited states for any 
N < 52 and R from 3 to 6 (for some selected N values we 
have also analyzed R values going up to N), using exact al- 
gorithms. For a given couple N, R we have determined all ex- 
act GS using a clever exhaustive enumeration scheme called 
b ranch- and-bound (b&b). In b&b, the problem is solved us- 
ing recursion. In a branching step, one of the variables oi is 
chosen, and it is eliminated by creating two subproblems in 
one of which <j l = +1, and in the other <Xj = — 1. The latter 
are solved recursively. A subproblem is solved by determining 
upper and lower bounds on its optimum solution value. The 
energy value E uh of the best known configurations serves as 
upper bound. The latter are updated whenever configurations 
with energy E ub or better could be determined. If a subprob- 
lem' s lower bound attains a higher value than E ub , no con- 
figuration with lower energy can be contained in it. Thus, it 
can be excluded from further consideration (fathoming step). 
For designing a practically effective algorithm, it is crucial to 
employ a strong lower bound with which subproblems can be 
fathomed early, keeping their total number reasonably small. 
For additionally determining all excited states with value at 
most x% away from optimality, we fathom a node only if its 
lower bound is worse than (1.0 + -^)E°, being E° the exact 
GS energy (determined in an earlier run of b&b). Mertens 11411 
has developed a b&b approach for determining exact GS for 
the model restricted to R = N. He could solve the problem 
up to N = 60, which marks the world record. Following Hill , 
we first narrow the search space by exploiting the fact that if 
the configuration a\ ... a n is a GS also its reversal er n . . . <j\ 
is a GS. The same is true when each odd spin and/or each even 
spin is multiplied by —1. We restrict ourselves to representa- 
tives of these symmetry classes. 

A lower bound on Hi is given by minimizing each 
C(d,i,R) 2 independently. We briefly sketch how we es- 
timate the value of such a minimum. In a subproblem, 
several spins might be fixed, all other spins are yet free. 
The already fixed parts give some contribution, say Z. We 
look at all partial sequences of free spins at distance d, i.e., 
Ci+fcd, (? i+ (k+i)d, Pi+id that are framed by fixed spins. 
Depending on whether Z is negative or not, the smallest pos- 
sible value of C(d,i, R) is either achieved when the free spins 
are all set equal or all alternating. The bound that we calcu- 
late from this is stronger than the one presented in [14] in the 
sense that we need to enumerate considerably less subprob- 
lems. For small and medium R we get better performance 
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if we start branching by fixing the spins in the middle of the 
sequence, expanding towards the boundaries, than if we start 
branching on the spins along the boundary, moving 'inwards'. 
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FIG. 1: GS energy as a function of R, for different N values. 
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FIG. 2: R=5. Average energy of states at one spin flip from the GS. 
The energy gap in Hi is 4. In the inset the probability distribution of 
these energies for TV = 50 (dashed) and TV = 51 (continuous). 

In Fig. [1] we show the GS energy as a function of R. Dif- 
ferent full lines correspond to different TV values (increasing 
from bottom to top). The lower dashed line joins data points 
with R = N. For TV 3> R (upper lines for low R values) 
the data accumulate on a limiting curve: the thermodynamic 
model with R fixed is well defined. For R < 12 we worked 
out an iterative procedure for computing the GS configura- 
tions for any TV value. The TV ;§> R limiting curve has, like 
the infinite range Bernasconi model, an erratic R dependence, 
but smoother than in that case. 

We have run extensive simulated annealing experiments for 
R < 13 (to be better discussed below) and we have identified 
models where the energy relaxes fast to the GS one, e.g. R = 
3, 6, and models showing a very slow relaxation, e.g. R = 
5, 7. In the rest of the Letter we focus on the R — 5 model, 
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FIG. 3: Probability distribution of the overlap among the GS and the 
first excited states. TV = 50 (left) and TV = 51 (right). 



which seems the best compromise between very short range 
interactions and glassiness on quite long time scales 01711 . 

In Fig. [2] we show the average energy difference S of states 
at one spin flip from the GS as a function of TV -1 . The limit 
N — > oo can be estimated very reliably, and is close to 2.8, 
i.e. more than ten times the gap (that in these units is equal 
to 0.2). In the inset we plot the probability distributions of 
these energies for TV = 50 and TV = 51: they make clear 
that the large mean value of 5 comes from a quasi totality of 
configurations that have energies that are much larger than the 
GS one. The distributions shown in the inset only depend, for 
large TV, on TV being odd (where there is a single GS) or even 
(where there is a TV/2 + 1-fold degeneracy of the GS). This 
is a strong evidence that GS are surrounded by high barrier of 
the order of several energy gaps. This property is shared also 
by low energy configurations, and it is what makes the energy 
landscape somehow golf-course-like. 

In Fig. |3] we show the probability of the overlap q = 

lEi cr i S ' (T !™ ST ' among the GS and the first excited 
states, for TV = 50 and TV = 51. We consider all first excited 
states and compute q with the GS that is closest in Hamming 
distance. The qualitative difference of the two distributions 
for q smaller than 1/2 is connected to the fact that the GS is 
not degenerate for odd values of TV. These probabilities are 
peaked at a value of q different from one, and by increasing 
TV the mean value of q stays well below 1. The number of first 
excited states very similar to the GS is small; they are typi- 
cally far from the GS. This is a further hint towards a glassy 
nature of the system. 

We discuss now the finite T properties of the model (mainly 
for R = 5). The system size is always TV = 10 6 . In the main 
panel of Fig.|4]we show the relative difference A between the 
system energy and the GS energy during very slow annealing 
and heating experiments (marked by leftward and rightward 
arrows respectively) [18]. Temperature values look small, but 
this depends on the factor 1/(R(R — 1)) used in Eq.©: the 
relevant energy scale is the gap, which is 0.2 for R = 5. While 
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FIG. 4: Cooling and heating experiments: for R = 6 the model 
fastly relaxes to the GS, while for R = 5 it has a glassy dynamics. 
The dashed line is the exact energy for R — 5. Inset: integrated 
response versus correlation (t w values as in Fig. [5}. 
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FIG. 5: Correlation function C(t m ,t w + t) as a function of t (inset) 
and of t/t w (main panel). Values of t w are 10, 10 2 , 10 3 , 10 4 , 10 5 . 



the R = 6 model fastly relaxes to the GS energy, the R = 5 
one shows a very slow relaxation: a tentative extrapolation to 
the adiabatic limit using an inverse power of the running time 
returns A ~ 0.048, i.e. is roughly 5% above the GS energy. 
Results from heating experiments look like a crystal melting, 
although the dependence on the heating rate is strong. For 
R . = 5 the exact energy is plotted with a dashed line: it is 
worth noticing the presence at T ~ 0.038 of a secondary peak 
in C v , that seems to enhance hysteresis in temperature cycles. 

We show in Fig. [5] the two-time autocorrelation function 
C(t w ,t w +t) = N^ 1 &i(t w )(Ti(t w +t) for different values 
of the waiting time t w = 10, 10 2 , 10 3 , 10 4 , 10 5 . When plotted 
as a function of t (inset), the aging behavior is clear, and very 
similar to the one observed, for example, in a Lennard- Jones 
mixture J2]. The oscillations are maybe due to the determin- 



istic nature of the model as in Ref. [7]. There is a strong de- 
pendence of the decay rate on the waiting time t w . We show 
in the main panel how data collapse when plotted versus t/t w : 
again, we have a very good agreement with the behavior ob- 
served in glassy systems. 

In order to test more quantitatively the out of equilibrium 
regime, we have also measured the integrated response to an 
infinitesimal field switched on at time t w . We have used the al- 
gorithm described in Ref. [15]. We show in the inset of Fig. @] 
the usual plot of the integrated response x(tw, t w + t) versus 
the autocorrelation C(t w ,t w + t) parametrically in t: in the 
region where C is not too close to 1 the data follow a line of 
slope smaller than 1 (in absolute value) as for the p-spin model 
and for glass-formers. 

We have defined a class of models that are potentially good 
descriptions of glasses. We have shown that already one of the 
simplest models of our class, the D = 1 and R = 5 model, 
has glassy properties. We have analyzed the low energy land- 
scape (introducing a new effective bound in the optimization 
process), and used Monte Carlo dynamics to qualify its finite 
T behavior. There is much interesting work left, on the mathe- 
matical analysis of the model, on the study of different values 
of R and the Kac limit, and on the D > 1 problem. 

We acknowledge interesting discussions with A. Billoire 
and S. Franz. We are aware that T. Sarlat, A. Billoire, G. 
Biroli and J.-R Bouchaud are studying a local version of the 
ROM model. 
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